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Abstract Orthogonal Fractional Factorial Designs and 
in particular Orthogonal Arrays are frequently used 
in many fields of application, including medicine, en¬ 
gineering and agriculture. In this paper we present a 
methodology and an algorithm to find an orthogonal 
array, of given size and strength, that satisfies the gen¬ 
eralized minimum aberration criterion. The methodol¬ 
ogy is based on the joint use of polynomial counting 
functions, complex coding of levels and algorithms for 
quadratic optimization and puts no restriction on the 
number of levels of each factor. 
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1 Introduction 


In this paper we present a methodology to find one of 
the best orthogonal arrays for the genera lized min imum 
a berra tion (GMA) criterion, as defined in Cheng and w] 
(|2004 h We refer to these designs as GMA-optimal de¬ 
signs. For an m-factor design, the GMA-criterion is to 
sequentially minimize the severity of aliasing between 
all the i-factor effects and the overall mean, starting 
from i = 1 (main effects) and finishing at i — m (m- 
factor interaction effects). 

The joint use of polynomial indicator functions and 
complex coding of levels provides a general theory for 
mixed le vel or thogonal fract ional factorial designs, see 
Pistone and Rogantin l2008ll. It also makes use of com¬ 


mutative algebra, see IPistone and Wvnn ( 19961) . and 


R. Fontana 

Department of Mathematical Sciences - Politecnico di Torino 
Tel.: +39-011-0907504 E-mail: roberto.fontana@polito.it 


generalizes t he appr oach to two-level designs as dis¬ 
cussed in Fontana et all (2000). This theory does not 
put any restriction either on the number of levels of 
each factor or on the orthogonality constraints. It fol¬ 
lows that our methodology can be applied to find any 
GMA-optimal mixed-level orthogonal array. 

Orthogonal Arrays (OAs) are frequently used in many 
fields of application, including medicine, engineering 
and agriculture. They offer a valuable tool for dealing 
with problems where there are many factors involved 
and each run is expensive. They also keep the statis¬ 
tical analysis of the data quite simple. The literature 
on the subject is extremely rich. A non-exhaustive list 
of references, mainly related to the theory of the de- 


of 

Bose (1947) and the following books: 

- r-r- 

Raktoe et al 

(1981), Collombier (1996), 

Dev and Mukerjee (1999), 

Wu and Hamada (2000). Mukeriee and Wu 20060 and 

Bailev (2008)). 


Orthogonal Arrays represent an important class of 
Orthogonal Fractional Fac t orial Desi gns (OFFDs), see, 
for example, Hedavat et all ( 19991) and Schoen et al ( 2010h . 
Indeed an Orthogonal Array of appropriate strength 
can be used to solve the wide range of problems related 
to the quantification of both the size of the main effects 
and the interactions up to a given order of interest. 

This paper is organized as follows: in Section [2] we 
briefly review the algebraic theory of OFFDs based on 
polynomial counting functions. The computation of the 
wordlength pattern of a design is described in Section [3] 
while we describe the algorithm in Section [4j Some ap¬ 
plications of the algorithm are presented in Section [5] 
Finally, concluding remarks are in Secti on E] Section [2] 
is closely based on Section 2 of Fontana (|2013l) . We in¬ 
clude it here to facilitate the reader. 
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2 Algebraic characterization of OFFDs 


In this Section, for ease in reference, we present some 
relevant results of the algebraic theory of OFFDs. The 
interested reader can find further inf ormati on, including 
the proofs of the pr oposit ions, in Fontana et all (200Q) 
and lPistone and Roeantin ( 2008h . 


2.1 Fractions of a full factorial design 


Let us consider an experiment which includes m factors 
Vj, j = 1,..., m. Let us code the n 3 levels of the factor 
T>j by the n^-th roots of the unity 


Vj — {o 


> ^n-j — 1 f i 


where w= exp (\f-f ^ k\, k = 0,..., nj — 1, j = 
1 ,..., m. 

The full factorial design with complex coding is V = 
Vi x • • • Vj ■ ■ ■ x V m . We denote its cardinality by ffV, 

#© = nr=i«i- 

Definition 1 A fraction T is a multiset (J 7 *, /*) whose 
underlying set of elements J 7 * is contained in V and /* 
is the multiplicity function /* : J 7 * —> N that for each 
element in J 7 * gives the number of times it belongs to 
the multiset T . 

We recall that the underlying set of elements J 7 * is the 
subset of V that contains all the elements of V that ap¬ 
pear in T at least once. We denote the number of ele¬ 
ments of a fraction T by ffF, with ffX = /*(0- 

Example 1 Let us consider m = 1, ni = 3. We get 

V = {1, exp y ^ , exp yj }■ 

The fraction T = {1,1, exp (ff—1 ^p)} is the multi¬ 
set (J 7 *,/*) where J 7 * = {1, exp (-v/^ pp)}, /*(1) = 2, 
and /*(exp (y/^1 ^ L )) = 1. We get #J" = /*(1) + 
/*(exp(^f)) =2 + 1 = 3. 

In order to use polynomials to represent all the func¬ 
tions defined over V , including multiplicity functions, 
we define 

— Xj, the j-th component function, which maps a 
point f = (Ci,..., £ m ) of V to its j-th. component, 

Xj : V 9 (Cl, • • • j Cm) 1 t C j G ■ 

The function Xj is called simple term or, by abuse 
of terminology, factor. 


- X a = •... • X“™, aeL = Z ni x---xZ nm i.e., 

the monomial function 


x a :v 9(Ci,...,Cm) -»■ cr-.-.-c* • 

The function X a is called interaction term 

We observe that { X a : a G L = Z ni x • • • x Z Um } 
is a basis of all the complex functions defined over V. 
We use this basis to represent the counting function of 
a fraction according to Definition [2] 


Definition 2 The counting function R of a fraction J- 
is a complex polynomial defined over V so that for each 
C G V, R( C) equals the number of appearances of C in 
the fraction. A 0 — 1 valued counting function is called 
an indicator function of a single replicate fraction J 7 . 
We denote by c Q the coefficients of the representation 
of R on V using the monomial basis {X a . a £ L}: 


R (0 = c « A '“(C). C G V, c a £ C . 

a eL 

With PropositionlTIfrom lPistone and Rogantin (2008), 
we link the orthogonality of two interaction terms with 
the coefficients of the polynomial representation of the 
counting function. We denote by x the complex conju¬ 
gate of the complex number x. 


Proposition 1 If J- is a fraction of a full factorial de¬ 
sign V, R = l c aX a is its counting function and 
[a—f}\ is the m-tuple made hy the componentwise differ¬ 
ence in the rings Z nj , ([au - ff } ni ,..., [a m - /3 m ]„ m ), 
then 

1. the coefficients c a are given by c a = ^" Q (C) ! 

2. the termX a is centered onT i.e., YIq&t Aa (C) = 
0 if, and only if, c a = C[_ a ] = 0; 

3. the terms X a and X@ are orthogonal on T if and 
only if, C[ a -p] = 0. 

We now define projectivity and, in particular, its re¬ 
lation with orthogonal arrays. Given I = {i \,..., ik} C 
{1,.. < ... < i k and C = (Ci, - ■ •, Cm) £ V we 

define the projection 7rj(£) as 

7T/(C) =0 = (Cii! • • ■ J C ik ) G Vi j x ... x Vi k . 

Definition 3 A fraction T factorially projects onto the 
/-factors, / = {ii,...,ik} C {1, i\ < ... < 
ik, if the projection ^/(J 7 ) is a multiple full factorial 
design, i.e., the multiset (2\ x ... x V ik , /*) where the 
multiplicity function /* is constantly equal to a positive 
integer over V ix x ... x V ik . 

Example 2 Let us consider m = 2, n\ = ri 2 = 2 and the 
fraction T = {(—1,1), (—1,1), 

(1,-1), (1,1)}. We obtain 7Ti(J 7 ) = {—1,-1,1,1} and 
^(J 7 ) = {—1,1,1,1}. It follows that T projects on the 
first factor and does not project on the second factor. 
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Definition 4 A fraction J- is a mixed orthogonal array 
of strength t if it factorially projects onto any /-factors 
with ffl = t. 

Proposition 2 A fraction factorially projects onto the 
/-factors, 

/ = {«i,..., 4} C {1,..., to}, h < ... < i k , if and only 
if, all the coefficients of the counting function involving 
the I-factors only are 0. 

Proposition [2] can be immediately stated for mixed 
orthogonal arrays. 


3.1 The wordlength patterm of a fraction 

Given a J 7 of the full factorial design T>, let us consider 
its counting function R = Yh a <=L c aX a . From item|T]of 
Proposition [I] the coefficients c a are given by 

C6.F 

or equivalently 

c « = ^E^o^(C)- 

w C ex> 


Proposition 3 A fraction is an orthogonal array of 
strength t if and only if, all the coefficients c a , a 
(0,..., 0) of the counting function up to the order t are 

0. 


3 Aberration criterion 


Using the polynomial counting function, Cheng and 33 
(2004) provide the following definition of the general¬ 
ized wordlength pattern ajr = (oq (J -),..., of a 

fraction T of the full factorial design V. 


Definition 5 The generalized wordlength pattern ajr = 
(ai(J r ),..., a m (iF)) of a fraction T of the full factorial 
design T> is defined as 



where ||a||o is the number of non-null elements of a. 

According to the algebraic methodology that we 
have described in Section 0 as c a ’s are complex num¬ 
bers, we should simply generalize Definition [5] as fol¬ 
lows. 


Definition 6 The generalized wordlength pattern ajr = 
(aq(.F),..., a m (J 7 )) of a fraction T of the full factorial 
design V is defined as 



where ||x ||2 is the norm of the complex number x. 

The generalized minimum aberration criterion is to 
sequentially minimize cni(J-) for i = 1 ,..., to . 

In Section we provide a formula to compute 
i = 1,..., to, given a fraction /CD. 


The square of the norm of a complex number x can 
be computed as xx. It follows that 

11 Ca 112 — Ca Ca 

To make the notation easier we make the non-restrictive 
hypothesis that both the runs f of the full factorial 
design T> and the multi-indexes of L = Z ni x • • • x Z„ m 
are considered in lexicographic order. We get 

(#©)IM2 = £*(0^ = 

C ec 

= (X^ T Y)(X^ T Y) = Y T X^X aT Y 

where X a is the column vector [£“ : £ S T>\, Y is the col¬ 
umn vector [R(() : ( £ V] and . T denote the transpose 
of a vector. We refer to Y as the counting vector of a 
fraction and we denote by H a = [ h t j : i,j = 1,..., #/?] 
the matrix X a X aT . By construction the matrix H a is 
Hermitian and positive-definite. 

Proposition 4 The square of the norm of c a is 

IM2 = j^y yTh r y 

where Z/j} = [Re(ft.ij) : i,j = 1,..., ffD] and R e(hij) 
it the I'eal part of the complex number hij. 

Proof For a quadratic form we have 

Y T H a Y = Y T (H a ) T Y 

The matrix H a is Hermitian: (Z/ a ) T = H a . It follows 
that 

m\K\\l = Y T H%Y 

where //)) = [Re(hij) : i,j = 1,..., ffD] and R e(hij) 
it the real part of the complex number hij. □ 

In this way we can compute the generalized word 
length pattern using only real valued vectors and matri¬ 
ces. In Proposition [5] we provide an explicit expression 
of the elements of the matrix Z/j}. 
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Proposition 5 The real part of the element hij of the 
matrix H a is 


cos 




i,3 = !)■•■,#£> 


From a computational point of view (see Section 
m it is useful to consider the Cholesky decomposition 
of the symmetric and positive definite matrix Hi 

Hi = Uj Ui to . 


where (zi,..., z m ) (resp. (t\,..., t m ) ) is the i-th (resp. 
j-th) row of L = Z ni x • • • x Z„ m and n is the lowest 
common multiple of n ±,... ,n m , n = lcm(ni,... ,n m ). 

Proof Let ( = (Ci,•••, Cm) be the i-th row of V. We 
have C k = exp(-y^l^- 2 fe), k = 1,..., m where (zi,..., z m ) 
is the i-th row of L = Z ni x • • • x Z nm . Analogously 
let r = (n,... ,r m ) be the j-th row of V. We have 
r fc = exp(</M ^tk),k = 1, ... ,to where (ii, ...,t m ) is 
the j-row of L. 

The complex conjugate of C k is Cfc = exp(— y/^l^-Zk) 
k = 1,..., to. It follows that hij can be written as 


Thus from Proposition [5] the wordlength pattern of 
a fraction P QT> can be written as 

a i(?) = \\ u * Y \\l i = 1,..., TO. 

3.2 GMA for mixed level orthogonal arrays 

From Proposition[3]we know that for orthogonal arrays 
of strength t all the coefficients c a of the counting func¬ 
tion up to order t (that is 0 < ||a||o < t) must be 0. 
The generalized wordlength pattern a_F of an orthogo¬ 
nal array P of strength t will be 


O-jp - 27T 

exp( v / —T— ai(t\— Zi)) ■... -exp(\/—T- a m (t m -z m )) 

Tli 


(0,..., 0, a t +i{P),a m (P)). 


or 


/ /— 2tt . n , , 

exp(v-l — (—ai(ti - zi) + ... + 

n ni 


~m ) ) ) 


It follows that the counting vectors must satisfy the 
following condition 

U,Y = 0 i = l,...,t 


where n is the lowest common multiple of ni,... ,n m , or, equivalently, 

n = lcm(ni,..., n m ). Taking the real part of hij we 

complete the proof. □ = ® 


Proposition 6 The generalized wordlength pattern ajr = 
(ati (P ),..., a m (P)) of a fraction P of the full factorial 
design D is 

= mry YT HiY i = 1,• ■ ■,m 


where H t = Y^\\ a \\ 0 =i H R- 
Proof From Definition [6] we have 



From item [Q of Proposition [Q we get Co = and 
therefore ||co||i = (#$) • We can also write 

(#P) 2 Y \\ Ca \\l = {#V) 2 Y Y t H%Y = 

||a||o=i IM|o=* 

= (#V) 2 Y t ( Y h r) y = 0 m 2 Y T HiY 

ll«l|o=t 

where Hi = Y\ a \\„=i H R. D 


where 


A* — 


C/i 

Ut 


_ Fontana ( 2013 ), Fontana and Pistonel ( 2013 ) and Fontana and Sar. 
( 2013h show that, given the full factorial design V , the 
counting vectors Y = [i?(C) : C, £ V] ol the orthogo¬ 


nal arrays P C V of strength t are the positive integer 
solutions of a system of linear equations, AY = 0. The 
present paper shows a different way to build the con¬ 
straint matrix A. 

Let us suppose that we are interested in orthogonal 
arrays of size N. The GMA-criterion will provide the 
OAs with the maximum strength t for the given size N. 
If the norm of U\Y is strictly positive, 11 C/iW11 2 > 0 for 
all the counting vectors Y it follows that no OA exists 
for that size N. 


4 An algorithm for GMA designs 

A counting vector of a GMA-optimal design can be ob¬ 
tained through the to steps of the algorithm below. 
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4.1 Input 

The input of the algorithm is made of: 

1. the number m of the factors and the number of level 
rii of the *-th factor, i = 1,..., to; 

2. the size N of the fraction T . 

4.2 Step 1 

Solve the following quadratic optimization problem 
min ||C/iT||! 

subject to . , 

' i t y = n, • (1) 

y = [yi\,yi e > 0 

Let us denote by Y* one solution and let Wf = 11 C/i lAj* 111. 


Table 1 Wordlength patterns of optimal fractions at differ- 


ent steps; m = 

5,711 = .. 

= n 5 = 

2, N = 16 


step 

ai (T) 

«2 (W) 

as (JO 

04 (T 7 ) 

05 OF) 

1 

0 

10 

0 

5 

0 

2 

0 

0 

1 

0 

0 

3 

0 

0 

0 

1 

0 

4 

0 

0 

0 

0 

1 

5 

0 

0 

0 

0 

1 


Table 2 Wordlength patterns of GMA-optimal fractions of 
different sizes; m = 5, 7 ii = . .. = 715 = 2, N = 16 


N 

ai(J0 

o 2 (T) 

as (JO 

a 4 (JO 

05 (J 7 ) 

Type 

6 

0 

1.11 

1.78 

1.44 

0 

OA( 6,2 5 ,1) 

8 

0 

0 

2 

1 

0 

OA(8,2 5 ,2) 

10 

0 

0.4 

0 

1.8 

0 

OA( 10, 2 5 ,1) 

12 

0 

0 

1.11 

0.56 

0 

OA(12, 2 5 , 2) 

14 

0 

0.2 

0 

1.08 

0 

OA(14, 2 5 ,1) 


4.3 Step k with k = 2,..., m 

Solve the following quadratic optimization problem 

'm.m\\Uk;Y\\l 
subject to 
1 t Y = N, 

< \\lhY\\l = W*, ■ (2) 

\\U k -iY\\ 2 2 = WU, 

y = [ yi\,yi e %,yi > 0 


5 Test cases 


We denote by OA(N, m-.. .-n m , t) a mixed level orthog¬ 
onal array with N rows, m columns (the i-tli column 
has n, levels, 7 = 1 ,..., m) and with strength t. 

The computations are made using 


one main module, written in SAS_IML, th at pre¬ 


pares the m optimization problems, S.AS Institute Inc 

( 120081 ) ; 


MOSE K that solves each optimization problem, MOSEK ApS 
(20141. 


Let us denote by Y£ one solution and let W£ = 11Vj.* 111. The simulation study has been conducted on a standard 
We observe that if W* = 0, j 6 {1,..., to} then the laptop (CPU Intel Core i7-2620M CPU 2.70 GHz 2.70 
condition W* = ||[/jU ||9 can be simply replaced by GHz, RAM 8 Gb). 

UjY = 0 . 


4.4 Output 

If \\UiY*h = ... = || 2 = 0 and \\U t+1 Y* || 2 > 0 

the solution Y r * of the last optimization problem, cor¬ 
responding to the m-th step, is the counting vector of 
an orthogonal array T C T> of size N and strength t 
that is optimal according to the GMA-criterion. The 
wordlength pattern of T is 

(0, • • •, 0, —^*11 w* + illl,..., j^y\\w* m \\D- 

If ||t/iF*|| 2 > 0 then the solution Y^ is a fraction that 
that is optimal according to the GMA-criterion but that 
is not an OA. 


5.1 Five factors with 2 levels each 

We consider five factors with 2 levels each. We initially 
choose N = 16. Using the algorithm described in Sec¬ 
tion [2j we obtain 5 fractions, that are the optimal so¬ 
lutions corresponding to steps 1,..., 5. The wordlength 
patterns of these fractions are shown in Table [TJ The 
GMA-optimal fraction, that is found at the last iter¬ 
ation, is an OA(16, 2 5 ,4) with a wordlength pattern 
equal to (0,0, 0,0,1). 

We now consider N = 6, 8,10,12. The wordlength 
pattern of the GMA-optimal solutions are presented in 
Tabled 
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Table 3 Wordlength patterns of optimal fractions; m = 


4, rii = 

2 ,n 2 = . 

..=714 = 

3, N = 

18 

step 

ai{T) 

o 2 (-F) 

03 (.T) 

04 (T) 

1 

0 

1.81 

1.09 

1.09 

2 

0 

0 

1.78 

0.22 

3 

0 

0 

0.5 

1.5 

4 

0 

0 

0.5 

1.5 


5.2 Four factors, one with 2 levels and three with 3 
levels each 

We consider four factors, the first with 2 levels and the 
remaining three with 3 levels each. We initially choose 
N = 18. We obtain 4 fractions, that are the optimal so¬ 
lutions obtained at step 1,... ,4. The wordlength pat¬ 
terns of these fractions are reported in Table [3] The 
GMA-optimal fraction, that is found at the last itera¬ 
tion, is an OA( 18,2 • 3 3 ,2) with a wordlength pattern 
equal to (0,0, 0.5,1.5). 


6 Conclusion 

The joint use of polynomial counting functions and 
quadratic optimization tools makes it possible to find 
GMA-optimal mixed-level orthogonal arrays of a given 
size. It is worth noting that the methodology does not 
put any restriction on the number of levels of each factor 
and so it can be applied to a very wide range of designs. 
The methodology works with the standard partition of 
the set of the monomial exponents, L = Z ni x Z„ m : 
main effects, 2-factor interactions, ..., m -factor interac¬ 
tion but it can also be easily adapted to work with any 
partition of L. The range of applications is limited only 
by the amount of computational effort required. 
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